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Abstract 

The structure of flame filaments resulting from chaotic mixing within a combustion reaction is con- 
sidered. The transverse profile of the filaments is investigated numerically and analytically based on a 
one- dimensional model that represents the effect of stirring as a convergent flow. The dependence of the 
steady solutions on the Damkohler number and Lewis number is treated in detail. It is found that, below 
a critical Damkohler number Da cr u, the flame is quenched by the flow. The quenching transition appears 
as a result of a saddle-node bifurcation where the stable steady filament solution collides with an unstable 
one. The shape of the steady solutions for the concentration and temperature profiles changes with the 
Lewis number and the value of Da cr u increases monotonically with the Lewis number. Properties of the 
solutions are studied analytically in the limit of large Damkohler number and for small and large Lewis 
number. 

PACS numbers: 82.33.Vx, 82.40.Ck, 05.45.-a 
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1 Introduction 



Many chemical and biological processes take place within an imperfectly mixed environment, examples in- 
clude the mixing of reactants within continuously fed or batch reactors [jO and in the spread of plankton 
blooms within oceanic currents |^|, [jj. In these situations the time-dependent fluid flow, within which the 
reactions are taking place, can lead to the chaotic transport of fluid elements Q, |], which become 
stretched into thin elongated filaments. This behaviour has been observed both experimentally by inserting 
dye droplets into the flow || |j| |l(], [n], [l^] and in numerical simulations Jll], . 

In a two-dimensional system, a convergent and a divergent direction can be assigned to any point in the 
flow associated with the eigenvectors corresponding to the negative and positive Lyapunov exponents —A 
and A of the chaotic advection [ll| . These directions are, respectively, tangent to the stable and unstable 
foliations of the advection dynamics. Any advected material line tends to align along the unstable foliation 
in forward time, or along the stable foliation in backward time. Thus, the stirring process smoothes out the 
concentration of any advected tracer along the stretching direction, whilst enhancing concentration gradients 
in the convergent direction, producing a quasi-one-dimensional 'lamellar' structure. 



This allows us to separate the original reaction-advection-diffusion problem along the (Lagrangian) 
stretching and converging directions. In the stretching direction any perturbation is spread by the advective 
transport. This is the dominant process in this case and is much faster than diffusion. In the convergent di- 
rection, however, all three processes of reaction, advection and diffusion are of equal importance and need to 
be considered together. Chemical reactions of the type A + B —* P taking place in one-dimensional lamellar 
systems were first studied numerically by Muzzio and Ottino Jl6|, [l?], |l^| and analytically by Sokolov and 
Blumen jl9). More recently, Clifford et al. investigated the evolution of a two-step (competitive-consecutive) 
chemical reaction in lamellar systems |^0|, . 

For a spatially smooth velocity field the relative motion of nearby fluid elements can be approximated by 
linearising the velocity field along the path of a fluid particle. Thus, in the convergent direction fluid elements 
approach each other with a relative velocity proportional to their separation, Sv ~ —a(t)Sr, where the rate 
of convergence generally fluctuates in time along the trajectory. The average behaviour can be approximated 
by replacing a(t) with its long time average given by A, the Lyapunov exponent of the chaotic advection. In 
the context of autocatalytic type reactions it has been shown, that the resulting one-dimensional 'Lagrangian 
filament model' can be used to describe the mean transverse profiles of filaments that propagate along the 
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divergent direction following the unstable foliation |2^, [23] . 

Fluid mixing has an important role in combustion [ p4| , |25| , |26| , p7| p8| . Our main aim is to consider the 
effect of chaotic mixing within a combustion reaction, which we model as a first-order process converting a 
fuel C to an inert product P through the reaction 

C^P rate = ck(T) (1) 

with exothermicity q. T is the (absolute) temperature and c is the concentration of reactant C . The 
temperature dependence of the reaction rate is given by an Arrhenius law with an ignition temperature T"j, 
namely 

where E and i? are the activation energy and the universal gas constant and ko is the (constant) pre- 
exponential factor. Our reason for choosing this form for the rate law is to avoid the difficulties associated 
with the cold boundary problem. An alternative approach is, rather artificially, to set the ambient tempera- 
ture to zero. The form given by (2) allows for a non-zero ambient temperature and the discontinuity in k(T) 
does not lead to any problems provided Tj is kept small. 

Here we are concerned with the structures of filaments that arise within the chaotic mixing and so, for 
present purposes, the nature of the time-dependent fluid flow is not important. This aspect will be described 
in a future paper |2{| along with details of the resulting chaotic transport and reaction processes. This leads 
us to consider the equations for the filament structure, 



r)T f)T\ rP'T 

P C p [--Xx-)= K ^ +q ck(T) (3) 



dc . dc d 2 i, 

di- Xx Q- x = D 8?~ ck{T) (4) 



on — oo < x < oo, t > 0. Here x measures distance transverse to the filament, t is time and p is the density, 
C p the specific heat, k and D are the thermal conductivity and diffusion coefficient, respectively. 

The convective terms in equations (3,4) can be interpreted as advection by a pure strain flow at a con- 
stant stretching rate —A (A > 0) along the convergent direction. We can associate A with the Lyapunov 
exponent of the chaotic advection. Equations similar to (3,4) were first proposed by Ranz l30[| and have been 
considered previously for somewhat simpler chemical systems [fl9| , fj"o| including autocatalytic reactions |p3) 
and excitable media p2| , where their significance in determining the nature of the chaotic mixing is brought 
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out. These previous studies considered only chemical systems with identical diffusivities for all species. An 
important new aspect investigated in this paper is the possibility of different diffusivities for heat and reac- 
tant, characteristic to the combustion process. 



We apply the boundary conditions 

T^T a , c^cq as|x|^oo (t > 0) (5) 

where T a is the ambient temperature (we are assuming that T a < Tj). For the time-dependent problem we 
initiate the reaction by applying a local temperature input (above Tj). 

To make equations (3,4) dimensionless, we introduce the variables 

— T — T a c _ , [X 

T= -, c=— , t = Xt, x = x\— (6) 

T b -Ta co V D y 1 

where T b is the burnt temperature 7& = T a + ^rr-. This leads to the equations, on dropping the bars for 



convenience, 



dT dT d 2 T 

t. x ~^~ = Le — ; 

at ox dx 



-x— = Le^- + DacK(T) (7) 



where 



^-x^ = ^-DacK(T) (8) 
dt dx dx 2 ' w 



Da = Le 



X ' pC p D 

are the Damkohlcr and Lewis numbers, respectively. The temperture dependence of the reaction (2) becomes 

Km = { cxp (r ^x=m+wj) > ifr>^ (9) 

1 0, if T < T 

where 

— ^ b R — — T~ — Tl ~ ^° 

£ ~ E ' ; T b ' ±l ~ T b -T a - 

The boundary conditions to be applied are that 

T^O, c->l, as|x|^oo (10) 

The Damkohlcr number Da characterises the ratio between the advective and the chemical time-scales. 
Large Da corresponds to slow stirring or equivalently fast chemical reaction and vice versa. This, together 
with the Lewis number Le, are our main bifurcation parameters. 



We start by considering the possible steady state solutions to equations (7,8) satisfying boundary condi- 
tions (10). 



4 



2 Steady states 



Here we consider the steady equations 

LeT" + xT' + DacK{T) = 

(11) 

c" + xc' - DacK(T) = 

where K{T) is given by (9) and where primes denote differentiation with respect to x. We look for symmetric 
solutions, applying the boundary conditions 

T'(0) = c'(0) = 0, T^O, c->l asx^oo. (12) 

When Le = 1 we can combine equations (11) to eliminate the reaction terms. Integrating the resulting 
equation and applying (12) gives T + c = 1 and then 

T" + xT' + Da(l-T)K(T) = (13) 

again subject to (12). 

Equations (11), or (13), were solved numerically using a standard NAG library routine for integrating 
boundary-value problems. The outer boundary conditions were applied at a large value of x allowing them 
to be satisfied with sufficient accuracy. This value depended on the width of the filament which, in turn, 
depended on both the Lewis and Damkohler numbers. The numerical integrations were found to be insen- 
sitive to point where the outer conditions were taken. The results are presented in figure 1, where we plot 
T profiles for Le = 1 (figure la) and both T and c profiles for Le = 0.1 (figure lb) and Le = 10 (figure lc) 
for representative values of Da (e = 1.0, (3 = 0.1, 7* = 0.001). In these figures the broken lines represent 
temporally unstable solutions (see below). As Da is increased for the stable solutions (full lines) the width of 
the filament increases with the temperature reaching a constant value within the central part of the filamant. 
This central temperature quickly reaches a value of unity for Le — 1.0, whereas, for Le = 10.0, it slowly 
approaches this value from below as Da increases. The reactant C is fully consumed in this central region 
of the filament. For Le — 0.1, the reaction zone is much thinner than for the other two cases. Much higher 
temperatures can be achieved within the filament, though these reduce to unity, now from above, as Da 
is increased. In all cases the unstable solutions (broken lines) reduce in both magnitude and extent as Da 
increases. 
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To see how the solutions change as a given parameter is varied, we need some measure of the solution. 
The profiles shown in figure 1 suggest that the integrated quantity 

/>oc 

It = / T(x)dx (14) 
Jo 

is a suitable measure. Note that, if we add equations (11), integrate the resulting equation and apply 

/•OO POO 

boundary conditions (12), we find that It = T{x)dx = (1 — c(x))dx. Graphs of It against Da are 

Jo Jo 

shown in figure 2 (for Le — 0.1, 1, 10). The figure shows that there is a critical value of the Damkohler 
number, Da cr i tl with a saddle-node bifurcation at Da cr a giving two solutions for Da > Da cr a and no 
solutions for Da < Da cr i t , apart from the trivial uniform solution, C{x) = l,T(x) = 0, corresponding to the 
extinction of the flame. Thus, the transition at Da cr u can be interpreted as a quenching transition when the 
burning of the fuel is not sufficiently fast to compensate for the diluting effect of the stirring, which tends to 
reduce the temperature of the filament. It increases with Da on the upper solution branch and decreases 
on the lower branch. This effect can also be seen in the profiles shown in figure 1. Note that the value of 
Da cr u appears to increase with Le. We also note, that the possibility of quenching a flame by a sufficiently 
strong shear flow has been shown recently in a rigorous mathematical study by Constantin et al [ p8[ . 

We can compute Da cr u as a function of the other parameters through looking for a non-trivial solu- 
tion to the homogeneous linear equations that are obtained by making a small perturbation to the solution 
of equations (11) (or (13)). This point will be made clearer when we discuss the stability of the steady 
states in the next section. A graph of Da cr n against Le computed in this way is shown in figure 3 (for 
e = 1.0, f3 = 0.1, Ti = 0.001). The figure shows that Da cr n approaches a finite value (approximately 8.5 for 
this case) as Le — * oo and approaches a non-zero value as Le — ► 0. 

The graphs shown in figures 1-3 show different forms of behaviour depending on whether Le is small or 
large. There is also a distinctive structure of the reaction zone when Da is large and it is this limit that we 
now consider. 

2.1 Solution for large Da 

Here we describe how the upper branch (stable) solutions behave for Da large. The form of the graphs 
shown in figure 1 suggest that, for Da 3> 1, the profiles have two regions of constant temperature and 
concentration. There is a central (fully reacted) region where T — T c , c = and an outer (unreacted region) 
where T = 0, c = 1. There then must be a reaction zone which joins these two regions and which smoothes 
out the discontinuities in temperature and concentration. To determine the structure of this reaction zone 
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we assume that it is centred oni = xo(Da). We then put 

X = Da}' 2 {x-x ) (15) 

and look for a solution by expanding 

T(X; Da) = T (X) + Da^T^X) + • • • , 
c(X; Da) = c (X) + Do _1 ci(X) + • • • , (16) 
T c (Da) = T c (0) + Da^T^ + ■■■ 

A consideration of the resulting equations suggests that x is of 0(Da 1 ^ 2 ). This suggests that we put 

x (Da)^Da 1/2 (a +a 1 Da- 1 + ■■■), Da > 1 (17) 

where the at are constants to be determined. 

The equations at leading order are then 

LeTg + a Tt i + coK(T ) = 

(18) 

c' ' + a c' - c K(T Q ) = 

subject to the matching conditions 

T — > 0, c — > 1 as X — > oo, r -»T c W, c ^0 as I^-co (19) 

Primes now denote differentiation with respect to X. If we add equations (18), integrate and apply the 
boundary conditions as X — > oo we obtain the equation 

LeT' a + a T + c + a c = ao (20) 

If we now let X — ► — oo in equation (20), we find, on using (19), that 

Tj°) = 1 (21) 

This result shows, from (6), that the central temperature within a filament is the burnt temperature Tb for 
large Damkohler numbers. 

It is the solution to equations (18a) and (20) subject to (19) and using (21) that determines a . This 
boundary-value problem was solved numerically in a way similar to that used for determining the general 
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filament solutions. A graph of a against Le is shown in figure 4 for e = 1.0, f3 = 0.1, Tj = 0.001. The figure 
shows that a , the position of the reaction front for large Da, increases with Le, in line with the profiles 
shown in figure 1. 

We now consider the equations at 0(Da^ 1 ). We find, as before, that we can add these equations to 
eliminate the reaction terms to obtain 

LeT'{ + a T{ + (X + ai)l£ + c'{ + a c[ + (X + ai )c' = (22) 

The boundary conditions to be applied are that 

7i -> T c (1) , ci -> as X -> -oo, 7i -> 0, Ci -> as X -> oo (23) 

If we integrate equation (22) and apply boundary conditions (19) and (23) we obtain, on using (21), 

/oo 
(1 - c - T )dX (24) 
-oo 

The integral in (24) can then be evaluated as (1 — Le)/ao by integrating equation (20). Thus we have 

T C W = 1 ~ 2 Le , with T c -l + Da- 1 ( 1 ~ 2 ) ■■■ as Da ( 25 ) 

Expression (25) gives — when Le = 1 and shows why, in the central region, the burnt temperature is 
approached from below for Le > 1 and from above for Le < 1. 

2.2 Solution for small Le 

Here we obtain a solution to equations (11) valid for Le <C 1. We assume that Da is of O(l) and start by 
making the transformation £ = Le _1 / 2 x. This leads to the equations 

T" + £T' + DacK(T) = 

(26) 

c" + Le (£c' — DacK(T)) = 

still subject to (12), where primes now denote differentiation with respect to £. Equations (26) suggest 
looking for a solution by expanding in powers of Le. However, we find that this is not sufficient and we need 
to expand in powers of Le 1 / 2 , namely 

T = T + Le 1 / 2 T 1 +LeT 2 + ---, c = c + Le 1 ' 2 ci + Le c 2 + ■ ■ ■ (27) 

At leading order c ' = from which it follows that 

c = 1 (28) 
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and then 

r o ' + tr o + DaK{T ) = (29) 
The numerical solution of equation (29), subject to (12), shows that it has a similar structure to the results 

poo 

shown in figures 1 and 2. From the co-ordinate transformation It = Le 1 / 2 It where It = 2o(£)^£- 

Jo 

A graph of It against Da is shown in figure 5. The main point to note about this figure is the existence 
of a critical Damkohler number Da cr i t — 3.0136 with two solution branches for Da > Da crit . This is the 
limiting value of Da crit at Le = in figure 3. Note that the value of It at Da cr i t decreases with Le with 
It ~ l.bmLe 1 / 2 as Le -» 0. 

At (^(Le 1 / 2 ) we have c" = 0. Hence C\ = £>i, where b\ is a constant to be determined. At 0(Le) we then 
have 

4 = DaK(T a ), 4(0) = (30) 
We can integrate equation (30) and using the result obtained from equation (29) that 

/>OC 

Da / K(T )d£ = It, we have that 
Jo 

c 2 J T £ + 5 2 as £ _> oo (31) 

for some further constant 6 2 . A consideration of the equations for the 7}, (j > 1) in expansion (27) shows 
that they can be solved so as to satisfy the outer bondary condition. 

Expression (31) shows that we require an outer region in which T = 0, in which x is now the independent 
variable and c = 1 + Le 1 / 2 c. At leading order we obtain 

c" + xc = 0, c ~ bi + It x + ■ ■ ■ as x — > 0, c^Oas noo (32) 

The required solution is 

exp(-s 2 /2)ds = &i 1 - erf(x/V2) (33) 



Using the matching condition for x small then gives b\ = — lyy ^- The concentration at the centre of the 
filament is then 

c(0) = 1- Le 1/2 I T (Da)^+---, as Le -» (34) 

We can now see why the terms of Le 1 / 2 are required in the expansion (27). Without these terms b\ would 
effectively be zero, the solution to equation (30) at 0(Le) would be the same and again would not satisfy 
the outer boundary condition ((31) still holds). An outer region would still be needed and be given by (32) 



though with bi — 0. A solution to this problem cannot be obtained along the lines given by (33). 



2.3 Solution for large Le 

The structure of the solution for Le ^ 1 is illustrated in figure 6 for Le = 10 3 and Da = 12.0 (e = 1.0, (3 = 
0.1, Ti = 0.001). The graph shows that there is an inner region around x ~ where c is very small and T is 
approximately a constant. We start our solution for Le large in this inner region. We leave x unsealed and 
look for a solution of equations (11) by expanding 

T(x;Le) = T c + Le- ( - 1+m '>T 1 (x)---, c(x; Le) = Le' m Cl (x) + ■ ■ • (35) 

where m (> 0) and T c 0) are to be determined. The equation for c\ is 

c'{ + xc[ - aa = 0, 4(0) = (36) 

where a — DaK{T c ) is a constant. Equation (36) can be solved in terms of Confluent Hypergeometric 
Functions j|l| as 

Cl {x) = AoexpH^)^^;!;^) (37) 
for some constant A . From (37) it follows that 

a;a ( 1 + '--) ( 38 ) 
The equation for T\(x) is T"(x) — —aa(x) which can be integrated to give 

Ti{x) = -^-^exp(-a ; 2 /2) 1 F 1 (^-; -; -) + B 1 (39) 

for some further constant B\. From (39) 

Ti(x) ; ^^XTs 7^" +2 (! + ■■■) as x^co (40) 

V ; 2(a + 2) (s±i)!2 a / 2 V ; V ; 

An outer region is required to satisfy the outer boundary conditions. The form of equations (11) suggests 

that the appropriate scaling for this outer region is £ = Le~ x l 2 x with c and T both of 0(1). Applying this 

ct 

scaling in expressions (38) and (40) shows that we should take m — —. With this, the problem for the outer 
region is, at leading order, 

T" + £T' + DaK(T) = 0, £c' - DaK{T) = (41) 

(primes denote differentiation with respect to £) subject to 

c ~ 3 £ Q + • • • , T ~ T c ^° r +2 + • • • as £^ 0, 

(a + 2)(a + 1) 

c -> 1, T -> as £ -> oo (42) 
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where A = 



2 a/2 ' 



The problem given by (41,42) has to be solved numerically and it is this solution that determines T c 

/>OC 

(and Aq). The results are shown in figure 7a, where we plot It = T(Qd£ against Da (for e = 1.0, (3 — 



0.1, Ti = 0.001). Note that It — Le 1 l 2 lT- The graph shows that there is a critical Damkohlcr number 



Da crit — 8.503, with two solutions for Da > Da cr i t - This value for Da crit is the limiting value for large Le 
seen in figure 3. The numerical solution gives T c and hence the exponent m = ex 
for expansion (35) in the inner region. Graphs of m and T c are shown in figure 7b. 



3 Stability of the steady states 



To determine the stability of the steady states, we make a perturbation of the form 



T{x,t) =T s {x) + e ult T 1 (x), c(x, t) = c,(x) + c wt ci(x 



(43) 



where T s and c s are the steady states (solutions of equations (11,12)) and T\ and c\ are small. We must 
allow for the possibility that uj and T\, c\ could be complex. Substituting (43) into equations (7,8) and 
linearizing gives 



subject to 



LeTi' + xT[ + Da( Cl K(T s ) + T lCs K'(T s )) - wTi 



c'{ + xc[ - Da( Cl K(T s ) + T lCs K'{T s )) - wci = 



T[(0) = ci(0) =0, Ti -> 0, ci -» as x -» oo 







(44) 



(45) 



Equations (44,45) are a homogeneous problem for the eigenvalue w and to force a non-trivial solution in the 
numerical integrations we also apply the condition Xi(0) = 1. 



The turning points in the It - Da curves (see figures 2,5,7) correspond to saddle-node bifurcations, i.e. to 
points where u) has a real zero. Thus solving equations (44,45) with u) = (and regarding Da as an unknown 
parameter) determines the value of Da at the saddle-node bifurcation, i.e. determines Da cr i t - This was the 
procedure that was used to determine the values of Da cr u shown in figure 3. We can use this procedure to 
examine how Da cr i t varies with the other parameters. We considered how Da cr u varied with e. To do so 
we needed to think about the way in which we defined our Damkohlcr number. Our form uses ko for the 
reaction rate, perhaps a more suitable form is to use the the reaction rate at the burnt gas temperature Xf, 
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to define a Damkohler number Da, related to our choice by Da — e _1 / c Da. This change is not particularly 
significant for the larger values of e that we used (e ~ 1) but can be significant for smaller values of e. In 
figure 8 we plot Da cr n against e for Le = 1.0, j3 = 0.1 (and Ti = 0.001). The graph shows that we can have 
nontrivial solutions even for small values of e though relatively large values of e are required to have reason- 
able values for Da cr it- For smaller values of e (e less than about 0.3) the values of Da cr it increase very rapidly. 

We solved equations (44,45) numerically for a range of values of Le, varying Da. In all cases we found 
that the lower solution branch was unstable, in fact on this branch w was real with ui > 0. The upper 
branch solutions were found to be stable in all the cases considered. Generally we found to to be real and 
negative, though for high values of Le (Le ~ 10 3 ) we found that u> could become complex though always 
had Re(uo) < 0. We were unable to locate any parameter values at which the solution became unstable 
(and oscillatory) through a Hopf bifurcation at Re(u) = 0. These conclusions were confirmed by numerical 
integrations of the initial-value problem (7,8,10). These showed that a threshold temperature input was 
necessary to generate a nontrivial solution and that, when this was applied, the corresponding upper branch 
solution was approached for t large when Da > Da cr it- For Da < Da cr u the numerical solution returned to 
the initial state and the reaction ceased as t increased for all inputs. 

4 Conclusion 

The main feature to emerge from our discussion is the existence of a critical Damkohler number Da cr u ■ An 
important aspect of our kinetic model (1,2) in the reaction-diffusion context is that any localized tempera- 
ture input has to be above some threshold value if spatial or statiotemporal structures are to be sustained 
at large times. The effect of the mixing by the flow, seen in our model (3,4) through the steady converging 
flow Xx, is to remove the heat from where it is generated by the exothermic reaction. This has the effect of 
increasing the threshold input needed for steady structures to form at large times. Thus, for small values of 
the Damkohler number, mixing is too strong, the heat generated by the reaction is dissipated, the reaction 
cannot be sustained and the system returns to its original state. This means that relatively large values 
of the Damkohler number are required for nontrivial states to form within the filaments. This appears at 
a critical value Da cr a of the Damkohler number where there is a jump from the totally quenched state 
(Da < Da cr it) to a fully developed state (Da > Da cr it). Thus the values of Da cr it found in the present 
study are an important guide in determining where sustained combustion will arise in the more general flow 
problem |||. 
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We found that Da cr i t is not particularly sensitive to the Lewis number (for the values of the other pa- 
rameters considered), see figure 3, increasing from about 3 if Le is put to zero and approaching a limit of 
about 8.5 for large values of Le. However, the resulting temperature and concentration profiles are strongly 
dependent on both Da and Le. For large values of Da a central, fully reacted core develops which is sepa- 
rated from the outer conditions by a relatively thin reaction zone. In this case combustion effects are spread 
out over a wide area of the filament. A similar situation arises for large values of the Lewis number. Here 
the influence of heat conduction is very much stronger than the diffusion of the reactant. The effect is to 
increase the spread of the reaction region to OiLe 1 ^ 2 ), with a thinner inner region of constant temperature 
T c somewhat less than the burnt gas temperature, see figures 6,7b. For small values of the Lewis number 
the reaction region becomes much thinner, large temperatures can be achieved (relative to the burnt gas 
temperature) and only small amounts of fuel are consumed. 

Criticality is an inherent part of the system and appears to be present for all values of the parameters 
associated with the Arrhenius kinetics. To get 'reasonable' values for Da cr u requires that e be relatively 
large (figure 8), which can be regarded as corresponding to smaller activation energies. This is in contrast 
to other combustion problems where criticality is seen only for the higher activation energies and disappears 
(sometimes through a hysteresis bifurcation) at a given activation energy. Finally we note that we found 
only stable steady structures when Da > Da cr i t ■ We found no further bifurcations to oscillatory behaviour. 

Finally, we note that the simple one-dimensional model (3,4) considered here is a somewhat idealised 
model of the full advection-reaction-diffusion problem in which the effects of fluctuations of the stretching 
rate, curvature of the filaments, interactions between filaments etc. are neglected. Such effects, however, 
may play an important role in certain situations and deserve further investigations. 
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Captions for figures 

Figure 1. (a) Temperature T profiles for Le = 1 (Da takes values from 6.965 to 100 for the stable so- 
lutions and from 7 to 75 for the unstable solutions). Temperature T and concentration c profiles for (b) 
Le = 0.1 (Da takes values from 5 to 70 for the stable solutions and from 6 to 20 for the unstable solu- 
tions) and (c) Le = 10 (Da takes values from 8.5 to 50 for the stable solutions and from 8.7 to 50 for the 
unstable solutions), (e = 1.0, (3 = 0.1, Tj = 0.001). The broken lines represent temporally unstable solutions. 

Figure 2. Graphs of I T (defined by (14)) against Da for Le = 0.1, 1, 10 (e = 1.0, (3 = 0.1, T l = 0.001). 
The lower branch solutions are temporally unstable. 

Figure 3. A plot of Da crlt against Le (e = 1.0, (3 = 0.1, T, = 0.001). 

Figure 4. The solution for Da >• 1; a , the thickness of the central, fully reacted region, plotted against 
Le, (e = 1.0, (3 = 0.1, T t = 0.001). 

Figure 5. The solution for Le <C 1; It, obtained from the numerical solution of equation (29), plotted 
against Da (e = 1.0, (3 = 0.1, T, = 0.001). 

Figure 6. The solution for large Le; concentration c and temperature T profiles for Le = 10 3 , Da = 12.0 
(e = 1.0, = 0.1, T l = 0.001). 

Figure 7. The solution for large Le; (a) It, (b) the central temperature T c (full line) and the exponent m 
(broken line) for the inner region expansion (35) plotted against Da, (e — 1.0, (3 = 0.1, Ti — 0.001). 

Figure 8. A plot of Da crlt against e (Da crlt = e" 1 / 6 Da crlt ) for Le = 1.0, (3 = 0.1 (T t = 0.001). 
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Figure 2: FiglbFucl 
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Figure 4: FiglcFuel 
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Figure 9: fig5 
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Figure 11: fig7a 
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